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The problem of the phonon bottleneck in the relaxation of two-level systems (spins) to a narrow 
group of resonant phonons via emission-absorption processes is investigated from the first principles. 
It is shown that the kinetic approach based on the Pauli master equation is invalid because of 
the narrow distribution of the phonons exchanging their energy with the spins. This results in a 
long-memory effect that can be best taken into account by introducing an additional dynamical 
variable corresponding to the nondiagonal matrix elements responsible for spin-phonon correlation. 
The resulting system of dynamical equations describes the phonon-bottleneck plateau in the spin 
excitation, as well as a gap in the spin-phonon spectrum for any finite concentration of spins. On 
the other hand, it does not accurately render the lineshape of emitted phonons and still needs 
improving. 
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I. INTRODUCTION 

The problem of the phonon bottleneck (PB) was rec- 
ognized in 1941 by Van Vleck who started the abstract 
of his article 1 with the sentence "The present paper is 
rather negative in character" . Analyzing different ways 
for the phonons emitted by magnetic impurities (hence- 
forth spins) to relax, Van Vleck found that in many typ- 
ical cases the phonon rates are insufficient to keep the 
phonon subsystem at equilibrium. The phonons emitted 
by the direct processes, forming a narrow resonant group 
in the energy space, should be absorbed by spins again, 
that drastically throttles the net relaxation rate of the 
latter. This happens, of course, if the concentration of 
spins is large enough. There is no PB for a single spin 
interacting with phonons in a macroscopic crystal, no 
matter what is the phonon relaxation rate. 

Van Vleck's "negative finding" excited other re- 
searchers for an activity that has not ebbed until now. 
Subsequent publications 2 ' 3 explored the analogy with 
trapping of resonant radiation in gases, considered ear- 
lier by Holsteini 4 - The corresponding theory deals with 
the spatial diffusion of emitted and reabsorbed photons 
or phonons from the body of the specimen towards its 
boundaries where they escape (see also the recent Ref. 
5) . Giordmaine and Nash 3 - cite a number of early exper- 
iments where the PB was observed. Recently indications 
of phonon bottleneck were seen in molecular magnets*^ 
The usual "fingerprint" of the PB is the decrease of the 
net spin relaxation rate if the thermal contact between 
the crystal and the holder is bad. Then the emitted 
phonons cannot efficiently escape from the crystal, and 
the only way of their relaxation are nonlinear phonon 
processes. The dependence of the PB on the thermal 
contact becomes pronounced at low temperatures, where 
the direct processes are dominating. At higher temper- 
atures the Raman processes strongly come into play. It 
was show n 8 ' 9 ' 10 that a kind of phonon bottleneck takes 
place for the Raman processes as well, but the effect is 
much weaker than that for the direct processes. 

Existing theories of the phonon bottleneck for the di- 



rect processes operate with diffusion, kinetic, or rate 
equations that are not derived from the first principles. 
These equations are set up using the balance of the ex- 
citation numbers and the energy between the spins and 
phonons, as well as the Fermi golden rule. The quasi- 
continuum of resonant phonons is considered in a sim- 
plified way as a single dynamical variable. One group 
of publications 3 -^ takes into account effects of spatial 
inhomogeneity in the sample. Another grou p 8 ' 11 : 12 ' 13 : 14 
further simplifies the problem by ignoring spatial effects 
and modelling the relaxation of phonons by a single rate 
of approaching the equilibrium at some bath tempera- 
ture. 

Related class of problems deals with the influence of 
spins on the phonon spectrum 1 ^ and with the resonance 
scattering of acoustic waves on spins (see Ref. [H and 
references therein). Jacobsen and Stevens 1 ^ shown that 
the hybridization between spins and phonons in a crystal 
with regularly spaced spins leads to a gap or stop-band in 
the phonon spectrum around the spin-phonon resonance. 
This can be considered as an extreme case of the phonon 
bottleneck, as, obviously, one cannot speak about a uni- 
directional energy transfer from spins to phonons, at least 
if the phonons are undamped. If spins are diluted, the 
resonance phonon gap should persist, in a reduced form. 
This effect was never discussed in the existing theories of 
the PB. 

Since Van Vleck's originating work, the main stress 
in the PB problem was put on the sufficiency or insuf- 
ficiency of the phonon relaxation processes to transport 
the energy away from the spins. A common idea is that 
the phonon bottleneck is something that happens if the 
phonons are not relaxing fast enough. One can find in 
the literature a dimensionless bottleneck factor that be- 
comes large if the phonon relaxation rate is small [see, 
e.g., Eq. (36) of Ref. Q. However, numerous publica- 
tions on the single- spin-plus-phonon-bath model do not 
care for the phonon relaxation at all, that is certainly 
correct. The problem arises only if the number of spins 
is macroscopic and their concentration is large enough. 
Then the PB occurs and, as a consequence, one has to 



2 



take into account the relaxation of phonons. This means 
that there should be another dimensionless bottleneck pa- 
rameter that is independent of the phonon relaxation or 
escape rate. 

The important component of the PB problem dealing 
exclusively with the energy transfer between the spins 
and resonant phonons has not recieved a due attention 
until now. Starting, for simplicity, with undamped har- 
monic phonons, one can ask what will be the accurate 
dynamics of the system, the large-time asymptote of the 
evolution, the energy distribution of emitted phonons 
under the bottleneck condition. One can expect that 
the system of equation s] 8 ' 11 ' 12 ' 13 ' 14 describing resonant 
phonons as a single dynamical variable follows from a 
more detailed energy-resolved theoretical framework via 
the integration over the phonon energies. 

Such a detailed description is likely to include an equa- 
tion of motion for the spin averages coupled to the sys- 
tem of kinetic equations for the phonons of all possible 
modes, considered separately. Are there any difficulties 
that prevent the derivation of such detailed spin-phonon 
equations from the first principles? It is puzzling why 
it has not been done yet, given that quantum kinetic 
equations result from the many-body quantum mechan- 
ics, under certain conditions. The procedure includes at 
first the derivation of the Pauli master equatio n 17 ' 18 that 
breaks the reversibility of the quantum mechanics. Then 
kinetic equations for particular observables can be ob- 
tained by averaging appropriate operators with the Pauli 
master equation. 

Studying the mechanism of the energy exchange be- 
tween the spins and the resonant group of phonons is 
the aim of this work. Identical spins S = 1/2 (two- level 
systems) randomly placed in the crystal are considered. 
One starts with undamped harmonic phonons and adds 
their relaxation at the very end in the same simplest way 
as was done in Refs. iEllllHll The main finding 
is that the narrow distribution of emitted phonons vio- 
lates the condition for the Pauli master equation to be 
valid. Hence no standard kinetic equations for the res- 
onant phonons can be derived. Instead, one obtains on 
this way equations with memory that are of little prac- 
tical use. In this situation it is better to step back and 
use a more basic system of dynamic equations including 
nondiagonal elements of the spin-phonon density matrix. 
Numerical solution of this system of equations shows a 
phonon bottleneck for the concentration of spins large 
enough. There is a nontrivial asymptotic plateau for the 
spin excitation, if the phonons are undamped. For a low 
level of the initial spin excitation the resulting system of 
equations can be solved analytically. The analytical solu- 
tion shows a resonance gap that is similar to that found 
by Jacobsen and Stevens 1 ^ for nondiluted spins. For the 
model with the damped phonons, the system eventually 
relaxes to the thermal equilibrium. However, the effec- 
tive relaxation rate for the spins is much less than the 
phonon relaxation rate in the case of a strong bottleneck, 
contrary to the first-glance expectation. 

The following part of the article is organized as fol- 
lows. Sec. [H] sets up the Hamiltonian of the spin-phonon 



system, its the quantum states and the Schrodinger 
equation. In Sec. IIIII the derivation of the Pauli mas- 
ter equation is reviewed and its applicability conditions 
are discussed. It is argued that for a narrow resonant 
group of phonons the standard kinetic formalism based 
on the master equation should fail because of the long- 
memory effect. Sec. IIVI contains the discussion of the 
conditions for the phonon bottleneck, and the dimen- 
sionless bottleneck parameter is introduced. Sec. [V] ex- 
plores the short-memory approach to the PB based on 
the master equation. One obtains a system of cou- 
pled kinetic equations for spins and phonons that con- 
tain powers of the energy ^-function and is thus incor- 
rect. It is shown how one can obtain previously pub- 
lished bottleneck equations 8 - ' 11 ' 12 ' 13 ' 14 by mathematically 
non-rigorous manipulations with ^-functions. In Sec. IVII 
one steps back to the master equation with memory and 
obtains a system of dynamical cqiations describing the 
emission/ absorption of phonons by spins and the bottle- 
neck. These equations are enchanced by including terms 
responsible for the relaxation of phonons in a simple way. 
Numerical solutions are presented and discussed. Further 
the analytical solution of the system of dynamical spin- 
phonon equations in the case of low spin excitation is pre- 
sented, and analytical results for the bottleneck plateau 
in the spin axcitation and the effective relaxation rate of 
the spin-phonon system under the bottleneck condition 
are obtained. Discussion is done throughout the article. 
Sec. I VIII contains a summary of the results obtained, as 
well as a discussion of the inhomogeneous broadening and 
the interplay between the phonon bottleneck and phonon 
superradiance* 1 ^ 

II. THE HAMILTONIAN AND SCHRODINGER 
EQUATION 

Consider a spin-phonon Hamiltonian for N$ two-level 
systems (spins) within an elastic body of TV cells 

H = H + V, (1) 

where 

Hp = — ^ a iz + hw-kola^ (2) 

i k 

describes spins and harmonic phonons, cr being the Pauli 
matrix. Neglecting the processes that do not conserve the 
energy, one can write this Hamiltonian in the rotating- 
wave approximation as 

v = E E K^ 01 4 + A^v) , 0) 

where = Vke _lk ri . Below V will be treated as a 
perturbation. The operator X 10 = cj_ brings the spin 
from the ground state |f) = |0) = (J) to the excited 
state ||) = |1) = (°) while A" 01 = a + does the opposite. 

To describe quantum states of the system of phonons 
and spins, we use a basis that is a direct product of the 
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states of all phonon modes and of the spins. The wave 
function can be written in the form 



*= E c {mjw}I1i^>IIk) 

{MW} i k 



(4) 



where fa = 0, 1 corresponding to the ground and excited 
states of the spins, respectively, and = 0,1,2,... are 
the occupation numbers of the phonon modes. We will 
use the shortcuts for different basis states 

{fa} = S, {pu}=V, {fa}{uw} = (SV) = W. (5) 

The Schrodinger equation (SE) for our spin-phonon sys- 
tem ifirff /dt = Hty can be written as a system of equa- 
tions for the coefficients cyy that can be obtained by act- 
ing on \1/ with the Hamiltonian of Eq. (JXJ) . The general 
form of the SE for the coefficients is 



■* d 

ih—c w 
at 
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VVi 



= e w c w +)^c Wl (W 
Wi 



V 



Wi 



where 



W H W) = hw 



(6) 



(7) 



is the unperturbed energy in the state W. It is convenient 
to introduce the slow amplitudes cyy via 



c w {t)e 



c w (t) 
then the SE becomes 



dt CW 



Wi 



W 



V 



(8) 



(9) 



Let us work out the concrete form of the SE, for the 
spin-phonon interaction V given by Eq. ([3]). In the ma- 
trix element ( W 



V Wij the state Wi differs from W by 

one spin flip and creation or annihilation of one phonon. 
Thus it is convenient to write Wi in the incremental form. 
In particular, Wi = W, lj, — lk means that in Wi the 
spin on the site i is excited and one phonon in the k- 
mode is annihilated, relative to W. This can only happen 
if /Ltj = in W. In the state Wi = W, 0j, +lk the spin on 
the site i is deexcited and one phonon in the k-mode is 
created, relative to W. This can only happen if fa = 1 in 
W. One obtains 



V 



r 



= 7=^kAiiV^iT 

v -/V 



1 



W,li,-lk) = --=^(1-^)^,(10) 



where fa and refer to the state W. With the help of 
this, Eq. ([9|) can be written in the form 



dew 
dt 



^ i k 

£ £ e^*-" ^ (1 - Mi) ^W,U,-i k (H) 



III. THE PAULI MASTER EQUATION 

Let us construct now the elements of the slow density- 
matrix for our closed system of spins and phonons by 
multiplying coefficients of the wave function 



Pw\W ~ c w c w 



(12) 



The equation of motion for them follows from Eq. © 
and reads 



d 

j t Pw\w 



£(w F Wi)e'( 

Wi 



-i[ui w i-LJ w! ]t_ 



lJ Pwiwi- ( 13 ) 



We restrict our consideration to the case of a nearly 
diagonal density matrix of the spin-phonon system. This 
will be the case if the spins are prepared in the initial 
state with random phases, so that the averages of the 
transverse spin components are zero, (S x ) = (S y ) = 0, 
and the phonons are at thermal equilibrium. We are go- 
ing to derive equations for the populations of the quan- 
tum states 



w 



Pw\w- 



(14) 



One can see from Eqs. (fTB")) with W' = W that diagonal 
elements Pw|w are dynamically coupled to nondiagonal 

elements such as Pwi|w- The latter are generated by V 
and thus are small. One can integrate them out using 
the equations similar to Eq. (Ti"3"|) . This yields 



PWi\WV-) — PWi\W (to) 



-Y 
-Y 



w 2 



W 



V 



V 



Wo 



dt'e l ( Ww ^ 



Jw ^'p W2 i W (t r ) 



-i(u> w -u) W2 )t~ 



(15) 



and /0 W | Wi (t) = fpyvi|>v(*)) ■ Remaining within the 

second order in the perturbation V, it is sufficient to 
drop the nondiagonal terms in the right-hand sides of 
the above equation. Plugging the resulting expressions 
for Pw 1 |vv(*) an d Pw\w 1 (t) into Eq. (fT3|) yields the mas- 
ter equation with memory 

^) = ^Ei 

Wi 

x cos [{u; w - wwj (* - t')} [P Wl (t') - Pw{t')} (16) 

In the usual case where the kinetic theory is appli- 
cable. 
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/ dt 








J t 



w 



V 



Wi 



and are smooth functions 

of Wi. The summation over Wi in Eq. (|16p goes over 
a wide energy interval limited by the maximal energy 
of the phonon reservoire fiw max in some terms and lim- 
ited by the thermal energy ksT in other terms. The 
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resulting expression is then peaked at t — t' < l/oj max 
or t ~ t' ^ h/{ksT). Both of these characteristic times 
are much shorter than the relaxation time 1 /Y for Pw (t) 
since the spin-phonon relaxation rate T oc |Vk| 2 is small. 
Thus one can make the short-memory approximation 
P w (t') => P w (t) and P Wl (t') P Wl (t) in Eq. ([IB]), 
After that the time integration can be done easily, 



*,,/ u \ ( i ,/m sin [(wwi - w w ) (* - *o)] 
at cos [(ww — %! ) (t — t Jj = 1 

t W Wl — W W 



In the kinetic time range 

lift 
£ - t ~ — > , 

r Wmax KflT 

one can replace 

sin[(a;w 1 - ui w ) (t - t )] 



(17) 



(18) 



7r5 (uw, — %) (19) 



that breaks the time reversibility of the quantum me- 
chanics and leads to the famous Pauli master equation 17 



— Pw = ^ ^wwi (Pwx — P\v) 

Wi 

with the Fermi-golden-rule detailed transition rate 
2tt / ~ v 2 



(20) 



r> 2 



w 



v 



Wi 



<5(wwi-w w )- (21) 



The sum of the latter over Wi 



Wi 



Wi 



(22) 

is the decay rate of the state W. One can write fyy as 



V" 



PW, (23) 



where iV is the total number of phonon modes in the 
system that is proportional to the number of atoms in it, 



w 



V 



Wi 



is the average over the resonant states 



Wi, and 



w) 



(24) 



is the density of states that satisfies J dujwp(ujw) = 1. 
The name "master equation" says that one can generate 
kinetic equations for different physical quantities from 
it by averaging appropriate operators over the quantum 
states W with Pyy. 

The reader can find a more extensive derivation and 
analysis of the Pauli master equation in Ref. [TH The au- 
thors argue that, as Eq. (f2TJ)) indeed suggests, a system 
initially in the quantum state VV will spread over all mu- 
tually accessible resonant states Wi in a nonoscillative 



way, the final result of the relaxation being the micro- 
canonical distribution Py^ 1 cc 6 (%, — %) . However 
plausible these arguments might appear, there is a prob- 
lem since one obtains a square of the energy (5-function 
in Eq. (|20|) in the case of the decay of an initial fully 
occupied state W: One (5-function is contained in the 
transition probability TyyvVi of Eq. (|21[) and another one 
is carried by Pyj 1 ■ In the derivation of the master equa- 
tion above [see discussion below Eq. (fi"6)) ] it was stressed 
that Pyvi( f ') should be a smooth functions of Wi for 
the short-memory approximation leading to Eq. (|20[) to 
be valid. Physically it means that the probability of the 
quantum system should be not sharp but distributed over 
many states with different energies. Only in this case one 
can rigirously derive kinetic equations. In the case of the 
spin-phonon relaxation via direct processes, the energy 
of the system is fixed and the emitted phonons build 
a narrow resonant group in the energy space that is as 
sharp as the energy (5-function describing the probability 
of this process. It will be demonstrated below that this 
results in the inapplicability of the kinetic approach to 
the description of the phonon bottleneck and that effec- 
tive bottleneck equations of Refs. l8l fli1fl2lfl3llT4l can only 
be obtained by mathematically incorrect manipulations 
with (5-functions. 



IV. THE ORIGIN OF THE PHONON 
BOTTLENECK 

As we have seen above, the Pauli master equation does 
not resolve lineshapes that is the origin of the difficulties 
of applying it to the PB problem. Finite linewidths, how- 
ever, follow from the time-energy uncertainty principle. 
In the simple case of one spin in a macroscopic crystal (or 
one atom in a free space or a large cavity) one can neglect 
the incoming term in Eq. (|20p . then the solution for the 
decay of the initially prepared state is Pw(t) — e~ Twt . 
The finite lifetime of the decaying state leads to the finite 
linewidth of the emitted phonons. At t — > oo from the 
Schrodingcr equation one obtain a 20 i 21 



Pw 1 - 



1 



v 



1 



(wwi — ww)' 



+ rfv/4 



(25) 



that satisfies J2 Wl Pyvi = 1- 

Apart of the Lorentzian dependence on the energy, 
PyVi also depends on the matrix element that can dif- 
fer for different directions of the phonon wave vector k. 
One can average Pyv x keeping the frequency a->yvi con- 
stant. Then using Eq. ([23]) one can rewrite the result in 
the form 



P Wi) — 



W 



/4 



nNp (uj w ) Y w ( %i _ Uw f + r^/4 



(26) 



that allows to estimate the probabilities of the states with 
emitted phonons. At resonance the Lorentzian factor 
here equals to one, so that the estimation is (PyVi) ~ 
1/iVr, where 



Nr = ttN p (%) Tyy 



(27) 
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The phonon bottleneck 




FIG. 1: The phonon bottleneck happens for B > 1 and it 
results in the reabsorption of the initially emitted phonons by 
spins. Fast phonon relaxation supresses the bottleneck. 



is the effective number of states in resonance with the 
state W, taking into account its finite linewidth. N-p ^> 1 
is the condition for using quasicontinuum approximation 
for the states of the system and replacing the sums over 
W by integrals. Only if Nr ^ 1 is fulfilled, one can write 
S {uj Wi - lu w ) in Eqs. (JUl) and 

For macroscopic bodies the factor N in Eq. (f27f makes 
Nr very large indeed. In this case {Pw} of Eq. is 
very small. This validates neglecting the incoming terms 
in the master equation that leads to Pw(t) = e~ Twt and 
Eq. ([25]) . The decay of an initially prepared state into 
the continuum can be understood from a very simple 
statistical argument. The excitation initially localized at 
the state Wo becomes equidistributed between the Nr + 
1 = Nr states as the result of the relaxation. Then the 
probability to find the system in any of the states 1/iVr 
is extremely small. Once the excited state of the spin 
has decayed, the excitation never comes back because it 
is statistically unprobable. 

For very small crystals the phonon modes can be so 
sparse that Nr is of order one or even smaller. If the 
spin is at resonance with only few phonon modes, the 
dynamics of the system is complicated and the spin can- 
not relax completely since it is being excited by phonons 
again. This is the case of the phonon bottleneck. The 
simplest case of the PB is that of the resonance between 
the spin and a single phonon mode that is described in 
quantum-mechanics textbooks as the resonance between 
two general quantum mechanical states. It is well known 
that the probability to stay in the intially prepared state 
is harmonically oscillating and there is no relaxation. 

The consideration above in this section pertains to the 
case of one two-level system (spin) relaxing to a phonon 
bath at zero temperature. Below we will consider, in 
particular, relaxation of N$ spins 1/2 prepared in their 
excited states and relaxing to a phonon bath at T = 
(see Fig. [T]). The initial state W (all spins excited, no 
phonons) decays, initially, into the states with one spin 
flipped and one phonon created. All arguments above are 
still valid for this process. However, the resulting states 
will further decay into states with more spins flipped and 
more phonons created. The number of emitted phonons 
can be by a factor N$ greater than that in the case of one 
spin. Then, in turn, these phonons will begin to excite 



the spin subsystem. In this situation one cannot neglect 
incoming terms in the master equation. A criterion for 
this can be obtained by the generalization of the statisti- 
cal argument above: For N$ spins at resonance with Nr 
phonon modes one can construct a dimensionless param- 
eter 



Nr irNp (wyy) Tw 

that controls the phonon bottleneck and can be called 
the bottleneck parameter. Note that here ns = N$ /N is 
the number of spins per unit cell, so that B is indepen- 
dent of the system size. If B < 1, then the most energy 
goes to phonons. If B > 1, then the energy mainly re- 
mains in the spin subsystem, if the resonant phonons do 
not transfer their energy elsewhere. The latter case cor- 
responds to the PB. It is clear that in most situations 
the inequality B ^> 1 should be fulfilled. Indeed, if one 
roughly replaces p (ww) =^ 1/^max, then B 3> 1 is equiv- 
alent to ryy/o-Wx "C ng. The left part of this inequality 
is very small, so that the inequality holds for any not too 
low concentrations of spins. 

It is tempting to take into account that the state result- 
ing from the decay of the spin, VVi , has its own relaxation 
rate Twi due to the damping of the emitted phonons. 
This would lead to the replacement Tyy Tyy + Twi in 
Eq. (|28p. where Tyy = T, the spin relaxation rate due 
to the phonon emission, and Twi = Tph, the phonon re- 
laxation rate. Then one can assume that the bottleneck 
should disappear in the case T p h 3> T. Below the reader 
will see that the sum T + T p h does not arise in the the- 
ory. The PB, indeed, is suppressed for large r p h, but the 
condition for this is more subtle and given by Eq. (|67|) . 

Inhomogeneous broadening that typically exceeds re- 
laxation rates tends to suppress the PB. Although it 
would be important for experiments, appropriate exten- 
sion of the theory will not be explicitly made here, for 
the sake of transparency. This extension is discussed in 
the concluding section of the paper. 

V. SHORT-MEMORY APPROACH TO THE 
PHONON BOTTLENECK 

In this section we try to derive equations describing 
the phonon bottleneck from the Pauli master equation 
Eq. (|20p , i.e., within the short-memory approximation. 
We will see that on this way one has to deal with powers 
of the energy 6- function that forces one to choose a better 
approach. 

A. The spin-phonon master equation 

Using the expressions for the matrix elements of Eq. 
(HfJ]), for the detailed transition probabilities of Eq. (|2~T|) 
one obtains 

27T 

fw;W,0i,+l k = -77 |Vk| 2 ^ (fk + l)S(u k - UJq) 
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fyV;W,l 4 -l k = ^ l^k| 2 (1 - Mi) ^k<5(^k - ^o) -(29) 

The master equation, Eq. (|20p. for the spin-phonon sys- 
tem then becomes 



X (Pw,0i,+l k - Pw) 



2tt 

at 



Ea-^Ei^i 2 ^^ 



X (f\ 



W,ii,-ik ~~ Pw 



(30) 



Let us now redefine so that they refer to the state 
in front of which they stand, instead of refferring to the 
state W, as initialy defined. Then the master equation 
takes the form 



dt 



2tt 



i k 

X [^k-Pw,0i,+lk - (^k + 1) Pw] 

EEl^^-^^K-^o) 



i k 



x {(vk + l)P W ,u,-i k -VkPw]- (31) 

This master equation can be simplified by taking into 
account that all N$ spins are identical and the state of 
the system depends on the global spin excitation number 



(32) 



rather than on the particular spin configuration {/^}. Re- 
calling the definition of the shortcuts for the quantum 
states in Eq. ([5]) one can replace Pw => P(sv) =^ ^m>C) 
and rewrite Eq. (f3Tj) as 



k 

X [^ k P (tl _ 1)(:Pi+llc) - (v k + 1)P^(V)] 
k 

x [(^ k + 1) P^+uv,-^) - WP,*,{V)] ■ (33) 



state. The average population of the phonon mode k is 
defined in a similar way as 



N s 

rik = E M! x _ /()! Z_ 

^=0 r v 15 ^' !/ k =0 q^kf q =0 



E^IIEV)' (35) 



It is never assumed that the spin and phonon parts of 
the density matrix factorize. 

It is convenient to introduce the distribution function 
for the spins 



N s l 



~ /-! 



M- (^s-M)! 



IIE^ 



(36) 



k i> k =0 



that is normalized as X)^=o = ■ Then Eq. (|34|) takes 



the form 



1 1 ^ 



(37) 



For a macroscopic number of spins Ns 3> 1 the distribu- 
tion function is sharply peaked at [i = /2. For instance, 
for uncorrelated spins one has 



N s l 



(38) 



wherefrom follows /i = /i max = pNs- Clearly / M remains 
sharp if there is a correlation between the spins via the 
emitted and absorbed phonons. 

One can also introduce the conditional probability n^.k 
as the total number of phonons in the k-mode in the spin 
state [i. It is defined by 



Nsl 



toKJH - ... 



Ii\ (N s - fi)l 



e^iiev). (39) 



so that, evidently, «k of Eq. (|35[) can be written as 

N s 

" k = E n M,k/M- ( 40 ) 

Since is sharply peaked, we will use 7ik = fi/i,k below. 
C. Kinetic equations for spins and phonons 



B. The observables 

The average excitation per spin is defined by 
1 1 Ns N I °° 

where the binomial factor accounts for the degeneracy of 
the microscopic spin states. One has p — if all spins are 
in the ground state and p = 1 if all spins are in the excited 



The time derivative of p defined by Eq. (|34|) can be 
calculated with the help of master equation Eq. (f3"3"|) . The 
right-hand side of the resulting equation can be simplified 
by introducing and n^k with the help of Eq. (|39|) . 
After some algebra one arrives at the equation 




+ ^Ei T/k i 2 ^ k -^)E^ 

b k p=0 
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x K*-i,k ( N S ~ M + 1) //i-i - ?V,kP/ M ] 

Ns 



2tt v 

NN S ^ 
6 k 



|Vk| 2 6 (wk - w ) ^ 



M=0 



x K+i,k (m+ 1)/, 



V,k (^5 - p) / M ] • (41) 



The first line of this equation does not contain the phonon 
occupation numbers and it describes the spontaneous 
emission of phonons. Here 

r = % E l^l 2 S (wk - w ) = 2tt (iVkl 2 )^ K) (42) 



is the single-spin decay rate, 

p(wo) = jr^Siu^-uo) 



(43) 



is the phonon density of states at the transition frequency 



(44) 



of the spins, and ( |Vk| / is the angular average 

dO k 



^1') = 



4tt 



Eq. (|4Tj) can be drastically simplified by shifting the /i 
index under the sum over p, so that only enters. As 
a result one obtains 



dp 
~db 



-rp + X! I y k| 2 $ (w k - wo) 

k 

X E /M ? V,k (iV s - 2/i) . 
ft=0 



(45) 



Here the sharpness of / M , see Eq. (|38| and the comment 
below, leads to the final simplification. In the sum over /i 
one can replace p =>• p = pNs and n^k =>■ nk, according 
to Eq. (|40|) , and then use the normalization condition for 
/„, see Eq. This leads to the final result 

§ = -r P + (i - 2 P ) ^ X! I yk l 2 5 ( Wk - w °) nk - ( 46 ) 

k 

The kinetic equation for the phonons can be derived 
in a similar way. The result reads 

^ = N S ^ |V k | 2 6 (wk - Wo ) [p - (1 - 2p)n k ] . (47) 

In fact, Eq. (|47|) could be guessed since, together with 
Eq. (UH]), it satisfies the excitation conservation 



pNs + ^ rik = const. 



(48) 



A disappointing feature or Eqs. f4"u]l and (|4"7|) is that pow- 
ers of <5 (ojk — Wo) enter both of them. If, say, there are no 
phonons in the initial state, then at short times rik will 
grow accordingly to Eq. |47|) as rik oc t5 (wk — u>o) ■ Then 
multiplication of this by S (lu^ — luq) leads to a mathe- 
matically incorrect expression. This problem was men- 
tioned already in the discussion of the applicability of 



the Pauli master equation at the end of Sec. The 
short-memory approach leading to the master equation 
is inapplicable for the description of the decay of an ini- 
tially prepared state, if the reabsorption processes have 
to be taken into account, as in the case of the phonon 
bottleneck. 



D. Ad hoc broadening of the 5-function 

Instead of stepping back to correct the error made in 
the derivation of Eqs. (|4l)|) and (|4"7| , one can choose a 
cheap solution using the regularization of the (5-functions 
by ascribing them a finite linewidth. One can, say, as- 
sume that the lineshapes are Lorentzian with a linewidth 
T, such as radiational decay lineshape given by Eq. (f2l)j) . 
In this case the square of the 5-function regularizes as 



2 1 

S (w k - w ) => — <5 (w k - w ) 

7Tl 



(49) 



Practically one can replace the energy (5-function in Eq. 
(HI as 



6 (w k - w ) =>• — , 

7Tl 



(50) 



because there is one more delta function in rik. In the 
simplest case of | Vk | 2 = \ V\ 2 independent of the direction 
of k and no phonons in the initial state one obtains 

| = _r p+ (i-2 P )^-L|vf £„ t 



-I>+(l-2p)— Ns(po-p), (51) 



where we used Eqs. (|42|) and (f27|) to transform 



2]Vf 
NT 



Nup (wo) 



r 



(52) 



and used the conservation of the excitation, Eq. (f48|) , po 
being the spin excitation in the initial state. Finally one 
obtains 



d r 

dt P = - Tp 



sr(p -p)(i-2p), 



(53) 



where the bottleneck parameter B is defined by Eq. 
This equation is a particular case of the bottleneck equa- 
tion that can be found in the Abragam & Bleany's book 13 
and it is similar to all other bottleneck equations & n ' 12 i 14 
published earlier and later. Non of these publications 
provides a derivation of the bottleneck equations so that 
it is impossible to judge whether these equations have 
been written ad hoc or derived in an incorrect way simi- 
lar to that described above. 
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VI. DYNAMICAL EQUATIONS FOR THE 
PHONON BOTTLENECK 



A. Derivation and analysis of the equations 

Having seen the origin of the breakdown of the kinetic 
description of the PB, one can easily correct the error by 
stepping back to the master equation with memory, Eq. 
([To]) . Then calculations similar to those of Sec. fVl yield 
the equations with memory for spins and phonons 

f = -|£l^| 2 ^Vcos[( Wk -c, )(i-O] 

x{p{t') + [2p{t')-l]n*{t')} (54) 



and 



<in k 



2Ws mi" 



N 



dt' cos [(wk — ujq) (t — t')] 



x{p(t') + [2p(t')-l}n k (t')}, (55) 



instead of Eqs. (|46|) and (f47|) . As there is no more the 
energy <5-function in these equations, one has to deal with 
all possible phonon modes. This means that in general 
the system of equations above should be solved numeri- 
cally. 

As equations with memory are inconvenient for nu- 
merical solution, it is better to remove the memory by 
introducing the additional dimensionless variable 



x{p(t') + [2p(t')-l}n k (t')}, (56) 



where Y is the spin-phonon decay rate given by Eq. (|42p . 
Evidently r k (£) is related to the nondiagonal element of 
the density matrix that was integrated out in Sec. IIIII 
One can see now that integrating out this nondiagonal 
element was unneccesary. The same results could be ob- 
tained directly from the Schodinger equation, Eq. (fTTjl . 
The resulting system of dynamical equations describing 
the direct spin-phonon processes has the form 



dp 
~dl 

dr k 
dt 
dn k 
~df 



TV 4-^ 



2|v k r 



Rer k 



(w k - w ) r k + T [p + (2p - 1) n k ] 



N s 2\V k \ 

n r 



Rer k . 



(57) 



Note that this system of equations is time-reversible, as 
the underlying SE, Eq. (fTTj) . However, Eqs. (f57|) are 
nonlinear since they provide a reduced description of a 
many-body quantum-mechanical system in terms of a few 
variables. We suppose that there is no spin-phonon cor- 
relation in the initial state and use the initial condition 



r k (0) = 0. 



(58) 



If | Vic | = |V| independently of the direction of the emit- 
ted phonons, then with the help of Eqs. (|42|) and (|28|) one 
can represent the coefficient in the third equation as 



N S 2\V\' 

n r 



YB. 



(59) 



This gives an idea of the strength of the PB in Eqs. (|57[) . 
In the case 6 « 1 and no phonons in the initial state, 
the generated phonon populations n k are small and can 
be neglected in the second of Eqs. ([57)) . After that the 
latter can be integrated and the result can be plugged 
into the first equation. Here one can make the short- 
memory approximation Re r k = rp7r<5 (o> k — luq) that is 
justified. This results into the well-known pure decay 
equation dp/ dt = —Yp. 

At this point one can include the phonon dissipation 
into Eqs. (|57[) in the si mplest pos sible way, generalizing 
the method of Refs. 0QI||lIlitl4|: 



dp 1 x - 

dt " 



dt 

dt 
dnk 
~df 



Rer k 



i (w k - u ) - -r ph 



ru + Y[p+ (2p- l)n k ] 



Y p h (n eq - n k ) 



N s 2 \V*\' 
N Y 



Rer k . 



(60) 



Here T p h is the phonon damping and n cq is the thermally 
equilibrium phonon population near the resonance, 



1 



exp 



\k B T ) 



1 



(61) 



At equilibrium one has n k = n cq , r k = 0, and p is ob- 
tained from the equation p + (2p — 1) n k = that yields 



1 



Poq 



1 + 2n e 



exp 



k B T 



1 



(62) 



an expected result. 

In the case P p h ^> Y the variables r k and n k can be 
expected to relax much faster than p. Thus after some 
time they should adiabatically adjust to the instanta- 
neous value of p. The same will happen in the case of 
strong bottleneck at asymptotically large times. Neglect- 
ing f k and n k in Eq. (|60|) one obtains 



n k 



|V k | 2 [p+(2p-l)n, 



eq J 



K-c r + ^-^iv k | 2 (2p-i) 



(63) 



for the adjusted value of n k . The second nonequilibrium 
term here is small and it has a Lorentz line shape with 
the phonon line width T p h if the last term in the denom- 
inator can be neglected. Note that for large r p h the dis- 
tribution of emitted phonons becomes a smooth enough 
function to make the short-memory Pauli master equa- 
tion applicable [see discussion below Eq. p^]) ]. In this 
case, however, the problem trivializes and the bottleneck 
disappears. One can see that the line width of emitted 
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FIG. 2: Numerical solution of the dynamical equations for the phonon bottleneck with nonrelaxing phonons at T — for 
different values of the bottleneck parameter B. (a) Low initial excitation of spins, po — 0.1. Asymptotes p(oo) given by Eq. (|74[) 
with no = are shown by dashed horizontal lines on the right, (b) High spin excitation, po = 1. In the case (b) the solution 
diverges for B > 100. The bottleneck plateau p(oo) grows with B but it should not exceed 1/2. 



phonons narrows in the case of the inverse spin popu- 
lation, p > 1/2. The reason is the stimulated phonon 
emission acting as phonon pumping that competes with 
their natural damping r p h. The strength of the effect de- 
pends on the bottleneck parameter B, as follows from 
Eq. ([M]). 

Similarly to Eq. (|63p . one can find the adjusted value 
Rer k and plug it into the first of Eq. JBDJ). If |Vk| 2 = \V\ 2 
independently of the direction of the emitted phonons, 
integration over k yields the equation 



J = -r»[p+(2p-l)n eq ] I 



where 



2B (1 - 2p) r 2 



(r P h + r; h (p))r; h (p) 



and 



r; h ( P ) = yJrl h + 2B(i-2 P )r*. 

One can see that the condition 

r ph > Vbt 



(64) 



(65) 



(66) 



(67) 



is needed to ensure T*(p) = T and thus to suppress the 
bottleneck. For B>1, this is a stronger condition than 
the first-glance expectation T p h > I\ If Eq. (|67|) is satis- 
fied, then the well-known solution of this equation 



P(t) = Pcq+(PO-Peq)e' 



-V T t 



r(l + 2n oq ) (68) 



is recovered. In the case of strong bottleneck, opposite 
to Eq. (f6"T)) . simplification of Eq. (|65[) yields 



P h 



y/2B(l-2 P ) 



(69) 



0.06- 



0.04 




400 t 500 



FIG. 3: Numerical solution for the spin excitation p(i) of the 
dynamical equations with relaxing phonons, Eq. (|60[1 , at T = 
for B = 10 and different values of the phonon relaxation rate 
r p h. Note that T p h S> T is needed to suppress the bottleneck, 
if B » 1. 



for p < 1/2. One can see that for strong bottleneck the 
effective spin relaxation rate does not depend on T at all. 
This means that the transfer of energy from the resonant 
phonons to the rest of the phonon bath or elsewhere is 
really the bottleneck of the whole process. On the other 
hand, the effective spin relaxation rate is not just r p h, as 
one could assume. 

It should be stressed that in the case B ^> 1 and 
Eq. (|67|) not satisfied, Eq. (f64|) becomes a reason- 
able approximation at asymptotically large times only. 
At small and intermediated times the full solution of 
Eq. (|6T)|) shows oscillations of p similar to probability 
oscillations in the soluition for two coupled quantum- 
mechanical states. The heuristic condition T p h > T 
mentioned at the beginning of this analysis is insufficient 
to establish adiabatic adjustment of rk and to the in- 
stantaneous value of p because of the large term ~ B in 
the third of Eq. (50]). 
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FIG. 4: Numerical solution for p(t) with relaxing phonons and 
knT / (huoo) = 2. Here the spins are initially in the ground 
state, po = 0, and they absorb the energy from phonons. 
For T p h = 0, depletion of the phonons prevents reaching the 
thermal equilibrium value of p (the dotted horizontal line) . 

B. Numerical solution of the dynamical 
spin-phonon equations 

It is not difficult to solve Eqs. (f5T|) and ([50]) numeri- 
cally. The results of the numerical solution of Eqs. ([57]) 
for macroscopic samples and no phonons in the initial 
state (i.e., T = 0) are shown in Fig. [J] The direct spin- 
phonon relaxation rate was taken to be T = 0.01 and the 
phonon cut-off frequency w max ~ 1, so that the short- 
memory approximation is valid in the absence of the bot- 
tleneck. For simplicity, the phonon modes have been dis- 
cretized equidistantly, ensuring Nr S> 1. Here and in all 
other numerical calculations we have set \V k \ = |V| 
for simplicity, to use Eq. (|59|) . One can see that for 
B = the solution for p(t) is the pure-decay expo- 
nential, p = poe~ rt practically in the whole time do- 
main. The initial quadratic dependence of pit) stemming 
from the time reversibility of Eq. (I57|) is confined to very 
short times of order l/w max and it is not seen in Fig. 
[2J The results in Fig. [2^ for the initial spin excitation 
Po = 0-1 show oscillating approaching a plateau that be- 
comes higher with increasing B. These oscillations is a 
memory effect that is absent in the earlier theories of 
the PBiMiii^ii^iii For low spin excitations, p« 1, the 
bottleneck equations can be solved analytically (see be- 
low) . The analytically obtained asymptotic values of p 
are shown in Fig. [5^ by dashed horizontal lines on the 
right. 

Numerical results for the initially fully excited spin sys- 
tem, po = 1, are represented in Fig. [2p. The results are 
qualitatively similar to those in the case of the low spin 
excitation, Fig. if the bottleneck parameter B is not 
too high. However, for B — 100 the solution of Eq. (|F7)) 
shows an unphysical divergence with p going to infin- 
ity. This indicates that our bottleneck equations are not 
accurate enough to describe the high-excitation strong- 
bottleneck case. High excitation means the inverted pop- 
ulation in the initial spin state, po > 1/2, that formally 
corresponds to the negative spin temperature. Asymp- 



totically the spin temperature should equilibrate with 
the temperature of the resonant phonon group. As the 
phonon energies are not bounded from above, phonons 
cannot have a negative temperature. Thus also spins can- 
not have a negative temperature in the asymptotic state. 
One can expect that for the initial spin inversion and 
huge values of B the asymptotic common spin-phonon 
temperature will be just very high, i.e., the asymptotic 
value of p will be close to 1/2. This means that the 
asymptotic value of p should saturate at 1 /2 with increas- 
ing B. This tendency is seen in Fig. [2p, only instead of 
the saturation at p — 1/2 one obtains a divergence. 

Figs. [3] and |4] show the numerical solution of Eqs. ([60]) 
for B = 10 and different values of the phonon relaxation 
rate T p h. In Fig. [3] the temperature of the phonon bath 
is T = and the starting spin excitation is po = 0.1. One 
can see that for all T p h > the spin excitation p relaxes 
to zero. The effective spin relaxation rate increases with 
r p h and asymptotically approaches T, as the PB is grad- 
ually suppressed by the phonon relaxation. However, the 
values of T p h neeeded to achieve the effective rate close 
to r are substantially greater than T, in accordance with 
the estimation in Eq. (|57)) . 

Fig. |4] shows the PB with the opposite direction of 
the relaxation, also with B = 10 for different values 
of r p h- Here po — and the temperature is nonzero, 
ksT / (Hu>o) = 2. Thus the spins absorb the energy from 
the phonon bath. If T p h = and B > 1, there are too 
little phonons in resonance with the spins, so that the res- 
onant phonons get depleted and the spins cannot reach 
the equilibrium value p eq given by Eq. (|62p . However, for 
all T p h > the phonons are replenished and the system 
reaches the equilibrium. Again, for T p h satisfying the 
strong inequality of Eq. (|67[) . the spin relaxation curve 
acquires its standard form without the PB, given by Eq. 

(EH). 

Studying the time evolution of the lineshape of emit- 
ted phonons, in the absence of the phonon relaxation, 
gives additional insights into the dynamics described by 
Eq. (|57[) . Fig. [5] shows the numerical results obtained 
with po = 1 and B = 0.1. The time dependence of n k 
at resonance in Fig. remains oscillating at very large 
times, where p(t) is already a constant. The lineshape 
rik at the maximal time of the calculation plotted vs 
the phonon-mode number (1000 corresponds to the reso- 
nance) is shown in Fig. [5p. It is strongly oscillating, and 
rik becomes negative for some modes at some moments of 
time. That is, the lineshape oscillates both as the func- 
tion of the phonon mode number and the time, as Fig. (5^ 
indicates. At large times, when the excitation of the spins 
has already reached its asymptotic value, the energy is 
still being exchanged between different phonon modes of 
the resonant group. Even for smaller values of B, the 
lineshape does not approach the radiational lineshape of 
Eq. (j2"5|) that is also shown in Fig. [5p. This implies that 
Eq. (|57| is still not accurate enough to correctly describe 
the lineshape of the emitted phonons. Still it yields some 
lineshape, whereas earlier theories ^ 11 ' 12 ! 13 ' 14 did not con- 
sider the lineshape at all. 
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FIG. 5: Emitted phonons in the absence of the phonon relaxation at T — for a moderate bottleneck, B — 0.1. (a) Time 
dependence n k (i) at resonance remains oscillating long after the spin excitation p(t) reached the bottleneck plateau; (b) 
Lineshape of the emitted phonons (at the maximal time) is narrow but not everywhere positive, that indicates that Eq. (I57|) 
needs improving. 



C. Bottleneck at low excitation 

In the case of low excitation of the spins throughout the 
process, pit) <C 1, Eq. (|6T))) can be linearized by setting 
p + (2p — 1) n k p — n k and then it can be solved by 
the Fourier transformation. We define the functions on 
the whole time interval (—00,00) and set them to zero 
at t < 0. Then one has to introduce the initial-condition 
terms PoS(t) and no5(t) into the first and third equations. 
Note that the initial phonon occupation can differ from 
the equilibrium value n eq . Fourier transforms are defined 
by 



f°° duj . 

P{t) = 2^ ( " )£ 



etc. With 



dtp{t)e lbJ \ 
(70) 

(71) 



k 

— = n Q 6(t)+Y ph [n cq 9(t)~n k } + — - g k 



Eq. ([60|) can be put into the form 



dt 

d^P 

~dt 



k - - (^k - W ) V>k - o r P h ^k + r b - ™k] 



k - (w k -w )£ k - ~r ph t/> k . 



(72) 



After the Fourier transformation one obtains a system of 
linear algebraic equations 

. . 1 ^ 2|F k |% 

- tup = Po - - ^ -^£ k 

k 



-wftk = n + r ph 



n. 



- 1 , N S 2\V^~ 
nk }> + — ^ £k 



-zw + 



at r 



-iw£ k = - (w k - uj q ) V> k - 2 r phCk + r [p - 
-iw^ k = (w k - w ) Ik ~ ^ r ph^ k - 



(73) 



Here one eliminates n k and ip k and expresses £ k through 
p. Then plugging £ k into the first equation yields an iso- 
lated equation for p that can be solved. In the resulting 
expression for p one sums over k. Making then the in- 
verse Fourier transformation one obtains, in particular, 
the asymptotic value 



p(oo) 



1 + VS 



(74) 



in the case T p h = 0, describing the phonon-bottleneck 
plateau. Also for r p h = one obtains 

c u\ t \ sin(ft k t) 
CkW = (Po - no) 7= — , (75) 



where 



( Wk -wo) 2 + sr 2 . 



(76) 



This result features the avoided level crossing between 
the energy levels of the spins and phonons, as £ k (t) that 
describes the dynamical spin-phonon correlations contin- 
ues to oscillate even at resonance, a-> k — loq. The splitting 
A = BY 2 in f2 k is proportional to the concentration of 
spins via B given by Eq. 1(28]) . One can see from Eq. (j5§|) 
that in the case of nondiluted spins, TVs = N, one has 



A 



BY 2 



2 V . This result is similar to that of Ref. 



For the model with relaxing phonons, in the strong- 
bottleneck limit B ^> 1 one obtains the following asymp- 
totic (r p h£ 3> 1) relaxation law 



= I Po + f^"o ) e 



-rt 



(77) 



12 



where 



r = 



rr 



2BT 2 



(78) 



is the effective relaxation rate of the spin-phonon system 
in the case of strong bottleneck, B 3> 1. One can see that 
for r ph < \f~BT the effective rate is f = Tph/v 7 ^ < 
r p h- The same result follows from Eq. (|69ll in the limit 
p -C 1. For <C 1 in Eq. (jUT)) the effective rate becomes 
r = T and the bottleneck disappears. 



D. Asymptotic stability of the bottleneck equations 

With the methods of the previous section it is pos- 
sible to analyze the asymptotic stability of the bot- 
tleneck equations with undamped phonons, Eq. ([57]) . 
Asymptotically p = const, so that one can substitute 
it into Eq. (|72p with r p h = 0. Assuming, for simplicity, 
| Vk | 2 = | | 2 = const and introducing a new variable k 



V^k = K + 



Tp 



Wk - w 



(79) 



one obtains the system of linear ordinary differential 
equations 



£k = - ( w k - w ) ^k + ( 2 i> - 1) r "k 

0k = ( w k - W ) Ck 
n k = BT^ k . 



(80) 



The eigenfrequencies of this system of equations are Slk 
and 



fi k = ±V(wk-w ) z + Sr 2 (l-2p), (81) 

c.f. Eq. ([75]) . For the inverse population, p > 1/2, the 
eigenfrequency f2k becomes imaginary in the frequency 
region near the resonance. This leads to an instability 
that is seen in the numerical results as positive divergence 
of p and thus negative divergence of rik near the reso- 
nance. The problem is that for the inversely-populated 
initial spin states, pa > 1/2 and B 1 the asymptotic 
value of p is close to 1/2, as discussed above. Then even 
small inaccuracies of the dynamic equations can render p 
slightly exceeding 1/2 that leads to the divergency. The 
fact that for po = 1 the instability happens only start- 
ing from B ~ 100 (see Fig. (2b) where p(oo) should be 
already very close to 1/2, indicates that the instability is 
driven by the inaccuracies of Eq. (|57[) . 



VII. DISCUSSION 

In the main part of this paper it was shown that the 
dynamics of two-level systems (spins) interacting with 
a continuum of resonant phonons via direct emisson and 



absorption processes, including the phonon-bottleneck ef- 
fect, is much more complicated than generally accepted. 
Since the emitted/absorbed phonon packet is narrow in 
the energy space, the usual kinetic description based on 
the short-memory Pauli master equation is invalid. Tak- 
ing into account memory amounts to stepping back to 
a dynamic description including nondiagonal elements of 
the density matrix. This cures the problem of the powers 
of the energy S- functions in the standard formalism while 
being capable of describing the PB. The resulting sys- 
tem of dynamical equations, Eq. (|57[) . can be enhanced 
to include the relaxation of phonons in a simple way, re- 
sulting in Eq. ([60]) . For the low spin excitation, p <C 1, 
these equations linearize and can be solved analytically. 
In a number of particular cases, such as small bottleneck 
parameter B of Eq. [)28p or fast phonon relaxation rate 
r p h, see Eq. (|67|) , simple results without a bottleneck are 
reproduced. 

The approach formulated here is not a final solution 
of the long-standing PB problem but rather a next step 
in improving the existing theories. A deficiency of the 
present approach is its inability to provide an accu- 
rate lineshape for emitted phonons, even the well-known 
Lorentzian shape of Eq. (|28p for the weak bottleneck, 
B <C 1. This should be the reason for the instability of 
Eq. (|57|) in the case of the strong bottleneck, B > 1, 
with spin inversion in the initial state. As the dynamical 
bottleneck equations have been obtained within the same 
method as the Pauli master equation, that is, by cutting 
the infinite chain of coupled equations for the nondiago- 
nal elements of the density matrix at the lowest possible 
level, an evident idea of improvement is to take into ac- 
count the next-generation nondiagonal matrix elements. 
This would result in a more complicated formalism, how- 
ever, that deserves to be dealt with in a separate work. 

In real systems the spin transition frequency luq is in- 
homogeneously broadened. That is, there are different 
transition frequencies u>i for different spins, distributed 
around the average value u>o with a width Aloq T. 
Obviously this leads to the reduction of the PB effect 
since more phonon modes can exchange excitation with 
the spin system, or, in other words, less spins can ex- 
change excitation with a given phonon mode. Thus one 
can expect that the situation will depend on the effective 
bottleneck parameter B ~ (T/Aujo)B <C B. 

Incorporating inhomogeneous broadening in the the- 
ory seems to be straightforward but it would result in a 
serious complication of the formalism. The problem is 
that spins cannot be described by a single variable p any 
longer. Spin with each particular cjj should exhibit their 
own dynamics, so that one has to search for the distri- 
bution p(u>i). It is not clear if any analytical solutions 
for the corresponding system of equations are available, 
and the numerical solution should be significantly more 
difficult than the one presented above. Note that inho- 
mogeneous broadening was mentioned in previous publi- 
cations on the phonon bottleneck (see, e.g., RefJ^) but it 
was treated as part of the "phenomenological" line width 
of the transition, without any frequency resolution of the 
type suggested in this paragraph. 
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Spin-spin interactions make the problem even more 
complicated. Although the effects of inhomogeneous 
broadening and spin-spin interaction are important for 
comparison with experiment, it does not make sense to 
consider them in the present paper for the reason stressed 
above. In principle, to be on a safe ground, one has to 
find a satisfactory solution for the core problem of the PB 
and obtain a correct line shape of the emitted phonons, 
before adding inhomogeneous broadening and other real- 
life effects to the model. 

Of course, one can argue that in the presence of inho- 
mogeneous broadening the phonon line shape should be 
dominated by the latter and one obtains a good-looking 
everywhere positive line shape by a kind of averaging over 
the spin transition frequency u>i, instead of the unsatis- 
factory line shape in Fig. [SJd. According to this logic, 
including inhomogeneous broadening would have more 
sense, as the next step, than struggling for the accurate 
line shape in the bare model. Practically this might be 
true, although one cannot easily find a microscopic jus- 



tification for such an approach. 

It should be stressed that the whole consideration 
in this paper pertains to magnetically diluted systems, 
in which the phases of phonons emitted by one spin 
and reaching another spin are random. This excludes 
magnetically dense systems such as molecular magnets. 
In the latter, emission of phonons by different spins is 
correlated, that in the ideal case leads to the phonon 
superradiance* 1 ^ The interplay of the phonon bottleneck 
and phonon superradiance is an exciting and challeng- 
ing problem, and its theoretical description requires in- 
cluding further correlators that have been ignored above. 
Obviously, for magnetically dense systems applicability 
of the ad hoc bottleneck equations ^ 11 ' 12 ! 13 ' 14 is even less 
justified than for the diluted systems. 
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